Drag Reduction by Microbubbles in Turbulent Flows: the Limit of Minute Bubbles 
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Drag reduction by microbubbles is a promising engineering method for improving ship perfor- 
mance. A fundamental theory of the phenomenon is lacking however, making actual design quite 
up-hazard. We offer here a theory of drag reduction by microbubbles in the limit of very small 
bubbles, when the effect of the bubbles is mainly to normalize the density and the viscosity of the 
carrier fluid. The theory culminates with a prediction of the degree of drag reduction given the 
concentration profile of the bubbles. Comparisons with experiments are discussed and the road 
ahead is sketched. 
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PACS numbers: 47.27.Rc, 47.55. Dz, 83.60.Yz 



The idea of reducing drag friction by placing a thin 
layer of air between a ship and its water boundary was 
patented already in the nineteenth century [lj. Drag re- 
duction by the injection of microbubbles into the tur- 
bulent boundary layer has been the subject of intensive 
research since the first experimental observation of this 
phenomenon in Q , and see the comprehensive review Q • 
The reduction of skin-friction drag by microbubbles has 
important technological and engineering advantages, es- 
pecially for marine transportation by huge and relatively 
slow ships like tankers, but also for many other appli- 
cations, such as hydro-foils, in-pipe transportation, etc. 
The voluminous literature on the engineering aspects of 
the problem cannot be referenced in full. It suffices to 
mention impressive results such as the microbubble drag 
reduction by about 80% on a flat plane, |j|, and up to 
32% on a 50 m long flat plane ship, see [j|. Some steps in 
understanding the phenomenon have been made. Ref |6| 
found that the drag reduction correlates with the maxi- 
mum void fraction in the boundary layer. It was under- 
stood that the "local distribution and shape [of the mi- 
crobubble void fraction C(r)] in the boundary layer have 
paramount influence in the drag reduction" pj . Many re- 
searches (see, e.g. (8J) found that effect of micro-bubbles 
decreases downstream and that the bubble size is another 
important factor influencing frictional resistance. 

Legner |9(] stated that the "decrease of the medium 
density as the gas concentration increases provides the 
primary drag reduction mechanism" . Unfortunately, the 
analysis of Ref. |9j does not contain any spatial depen- 
dencies, taking the distribution of bubble void fraction to 
be homogeneous. In addition Legner [9| concluded that 
the increase of the dynamic fluid viscosity, caused by the 
bubbles, leads to increase of frictional drag. In contra- 
diction, other studies (see, e.g. Ref. [lfJ() lead to the 
opposite conclusion: that the increase of the viscosity, 
caused by microbubbles, decreases the friction drag. To 
date this confusion has not been resolved theoretically. 

The aim of this Letter is to offer a theory of drag reduc- 
tion by microbubbles in the limit that their diameter d is 
very small (d — ► 0), and the void fraction C(r) is fixed, 
and not too large (C(r) < 0.1). In addition we will as- 
sume that the scale of variation £q = C(-r)/| VC(r)| <C z 
where z is the distance from the wall. The advantage of 



a theory in this limit is that we can show quite rigorously 
that the only mechanism for drag reduction available in 
this limit is provided by the reduction of the fluid den- 
sity and the increase in the fluid viscosity. This is not to 
say that there are no additional possible mechanisms of 
drag reduction by microbubbles due to their influence on 
the structure of turbulence, including near wall coherent 
structures [TI], [13, UM ■ The theoretical description of 
such effects is however very difficult; they stem entirely 
from finite bubble-size effects, and they should be taken 
only as a further step in the development of the theory. 
As a starting point for the theoretical development 
we take the two-fluid description of turbulent flows with 
bubbles which is presented in Ref. [bj|. In this descrip- 
tion the bubbles are of diameter d which is very small. 
We do not consider individual bubbles, but rather de- 
scribe them by a field of void fraction C(r,t) <C 1 and 
velocity w(r,t). The carrier fluid has density po, viscos- 
ity no and velocity U(r,t). We will take the air density 
of the bubbles to be zero and the acceleration due to 
gravity, g, to act in the z direction which is normal to 
the wall. Disregarding terms of the order of d 2 one write 
the equation of motion 



<'-*f 



18C> 
d 2 



(w - U) - (1 - C) Vp (1) 



+(1 - C) Po g + 2(1 - C)V • (jiEm) - -voCV 2 U 



2^o V ■ Ef-Vp 



^ (w-U) + \^ 2 U = . (2) 



In these equations 
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The effective viscosity which appears in these equations 
is determined by the bubble concentration, 



H= (1 + 5(7/2)^0 



(3) 



These equation should be supplemented with the conti- 
nuity equations 



d(l~C)/dt + V ■[(l-C)U] = 0, 



(4) 



dC/dt + V -(Cw) = 



(5) 



We now simplify the equations further in the limit d — > 
by evaluating the term proportional to d~ 2 in Eq. (JTJ 
using the same term in Eq. ^. We find 

(l-C)po^- = (l-C)[-V P + Po g + 2V-(pE m )} 
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In the same limit U m = U and E n 

After some further simplifications in which we retain 
only terms linear in C(r) one gets: 
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V • U = , DC/Dt = . (8) 

where the effective density of the suspension is 



p=(l-C)p 



(9) 



The important conclusion is that dilute (C <C 1) solu- 
tions of minute microbubbles (d — > 0) can be described 
by a one-fluid model with modified density p and viscos- 
ity p. Note that velocity field remains incompressible; 
this result is valid for minute microbubbles d — > for ar- 
bitrary concentrations C . Having these results at hand 
we are poised to offer a theory of drag reduction that is 
quite similar to the theory by the same authors for drag 
reduction by flexible polymers [15J. 

Consider a flow in channel geometry (with half channel 
width L); the mean flow is in the x direction, the wall 
normal direction is z and the span-wise direction is y. 
We take the bubble concentration C(r) to be given and 
time independent. The fluid velocity U(r) is a sum of its 
average (over time) and a fluctuating part: 

U(r,t) = V(z) + u(r,t), V(z) = (U(r,t)) . (10) 

For channel flows all the averages, and in particular 
V(z) => V(z), are functions of z only. The objects that 
enter the theory are the mean shear S(z), the Reynolds 
stress W(z) and the kinetic energy K(z); these are de- 
fined respectively as 

S(z) = ^-, W(z) = -p(z)(u x u z ), K(z) = ^(\u\ 2 ) . 

Under the assumption 4 «j/ we derive point-wise bal- 
ance equation for the flux of mechanical momentum, re- 
lating these objects [l5lllfil|. Near the wall it reads: 

p{z)S{z) + W{z) = p'L, forz<L. (11) 

On the RHS of this equation we see the production of 
momentum due to the pressure gradient; on the LHS we 
have the Reynolds stress and the viscous contribution to 
the momentum flux, with the latter being usually negligi- 
ble (in Newtonian turbulence p — p ) everywhere except 
in the viscous boundary layer. 



A second relation between S(z), W(z) and K(z) is 
obtained from the energy balance. The energy is cre- 
ated by the large scale motions at a rate of W(z)S(z). 
It is cascaded down the scales by a flux of energy, and 
is finally dissipated at a rate e, where e = p(z)(\Wu\ 2 ). 
We cannot calculate e exactly, but we can estimate it 
rather well at a point z away from the wall. When 
viscous effects are dominant, this term is estimated 
as [p(z)/p(z)](a/z) 2 K(z) (the velocity is then rather 
smooth, the gradient exists and can be estimated by the 
typical velocity at z over the distance from the wall). 
Here a is a constant of the order of unity. When the 
Reynolds number is large, the viscous dissipation is the 
same as the turbulent energy flux down the scales, which 
can be estimated as K(z)/t(z) where t(z) is the typi- 
cal eddy turn over time at z. The latter is estimated as 
yj p(z)z/by/K(z) where b ~ 1 is another constant. We 
can thus write the energy balance equation at point z as 
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K(z)=W(z)S(z), (12) 



where the bigger of the two terms on the LHS should 
win. We note that contrary to Eq. i|ll|) which is exact, 
Ea. (|12|) is not exact. It was shown however to give excel- 
lent order of magnitude estimates as far as drag reduction 
is concerned [la. Il7| . Finally, we quote the experimental 
fact |18L iigj l that outside the viscous boundary layer 



W(z) = c 2 K(z) 



(13) 



with the coefficient c rigorously bounded from above by 
unity (The proof is |c| 2 = \W\/K < 2\(u x u z )\/(u 2 x +u 2 z ) < 
1, because (u x ± u z ) 2 > 0). 

We can change variables now in favor of wall units 
according to 



S+ = p(z)S/p'L, z+ = zyJp{z)p'L/p{z) 
W + = W/p'L , K+ = K/p'L . 

In these units our balance equations read 
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This set of equations is readily solved, giving 

S + (z+) = l, forz+<z+, 



(16) 



S+(z+) = 
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In these equations we defined K = c 3 /b, z+ = a/c. The 
mean velocity anywhere in the channel can be obtained 
by integrating, 
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-S + (z+(z'))dz' . (18) 
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FIG. 1: Plots of the model response function \ and its contri- 
butions Xp an d Xa» due to the density and viscosity variations. 
These functions are computed from our model shear function 
S + (z + ) as shown in Eqs. JTHHIZJ plugged into Eqs. (EU E2J • 



A measure of drag reduction is the relative increase in the 
mean centerline velocity in the bubbly flow with respect 
to the neat Newtonian fluid: 



AV = V huh (L)-V N {L) 



(19) 



Clearly, AV > corresponds to the drag reduction, while 
AV < to the drag enhancement. We obtain an expres- 
sion for AV from Eq. (|18|l by expanding to linear order 
in C(z) (where our equations are valid anyway): 



AV+ = AV^/po/p'L = \ 
Jo 



x(z + )C(z+)dz+ . (20) 



Here the response function x consist of two parts, one 
due to the density variation Xp an d the other due to the 
viscosity variation x^'- 
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(21) 
(22) 

(23) 



In writing Eq. 1|20[1 we have used the fact that in exper- 
iments the bubbles tend to be localized in a finite region 
near the wall, i.e. C(z) — > sufficiently fast as z — > oo, 
and we extended the integration range to infinity. 

Eqs. (|20[l - i|23[l are the main theoretical predictions of 
this Letter. To complement the theory we present now 
estimates of the numerical value of the expected drag 
reduction, and compare it with a relevant experiment. 
The simplest model takes the parameters in Eq. (??) 
as 2-independent, and in agreement with the classical 
von Karman boundary layer theory, i.e. k w 0.436 and 
Zy ~ 5.6. Evaluating the response function x with these 
parameters results with the findings presented in Fig. ^ 



We see that in the viscous layer, where Eq. lfTB|) is rele- 
vant, Xp ~ while Xfi is negative. This means that hav- 
ing a bubble concentration in this region does not buy us 
drag reduction due to the density variation, but it leads 
to drag enhancement due to the viscosity increase. This 
is far from being surprising, since in this region the mo- 
mentum flux is dominated by the viscous term /iS. For a 
fixed momentum flux any increase in viscosity must de- 
crease S and correspondingly lead to drag enhancement. 
The most efficient drag reduction can be obtained by 
placing the bubble concentration out of the viscous layer, 
but not too far from the wall, say at 6 < z + < 30. In 
this region both the decrease in density and the increase 
in viscosity lead to drag reduction. The momentum in 
this region is transported mainly by the Reynolds stress 
—p(u x u z ). The effect of density reduction is absolutely 
clear: it leads to the reduction in momentum flux. For 
a given momentum generation p' L this has to result in 
the increase of the mean momentum of the flow. More 
interesting and counter-intuitive is the effect of increas- 
ing viscosity. In order to understand it, we remind the 
reader that for intermediate values of z + there is no well- 
developed turbulent cascade, and outer and inner scales 
of turbulence are of the same order of magnitude. There- 
fore the increase of viscosity reduces the turbulent energy, 
in contrast to fully developed turbulence where changes 
of viscosity simply modify the Kolmogorov scale without 
any effect on the turbulent energy, that is dominated by 
the outer scale. The decrease in turbulent energy here re- 
duces the Reynolds stress, see Eq. <|13fl . It is interesting 
to note that this effect of increasing viscosity is essentially 
the same as the mechanism for drag reduction in the case 
of elastic polymers [111 113 • For polymers, however, the 
increase in viscosity can be very significant and the linear 
approximation that is used here is not applicable. 

In comparing with experiments we need to consider 
low bubble concentrations. An interesting experiment 
was reported in [lj|, where both C(z) and the V + (z + ) 
are shown. We note that this experiment deals with a 
developing boundary layer rather than a steady chan- 
nel geometry, but near the wall the Reynolds number 
can be considered rather time- independent. Digitizing 
the published profiles C(z) and integrating them numer- 
ically against our function x( z ) we obtain results for AV^ 
which appear in good agreement with the data of |l0|. 
as long as C(z) is small, C < 0.1. For the two lowest 
values of bubble concentration we agree with the data to 
within 10-20%, which is definitely within the experimen- 
tal error bars. For higher values of C(z) the results of 
the experiment become sensitive to nonlinear effects. 

It appears extremely worthwhile to test the theory pre- 
sented here by numerical simulations that would be de- 
signed to do so. We should stress that a careful measure- 
ment of S + (z + ) in cither experiments or simulations, in 
addition to a determination of AV, can provide a very 
good test of our theory. Eq. I|20|l is more general than 
our model i|17|) . and it can be tested directly if S + (z + ) 
and its z + derivative are known. Since the response func- 
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FIG. 2: Plots of the response function \ and its contribu- 
tions Xp an d Xp due to the density and viscosity variations 
computed from the simulation data of 20]). 



dictions of our model for x are excellently reproduced by 
the numerical data, even though the smoother cross over 
between the viscous and logarithmic layers translates to 
smoother functions Xp an d Xm ■ A similar comparison for 
channel flow with bubbles will shed important additional 
light on our approach. 

We reiterate that additional nonlinear contributions to 
drag reduction are expected to come in when the concen- 
tration increases, and especially when the bubble diam- 
eter d increases. One should definitely examine theoreti- 
cally the nonlinear and finite size effects and incorporate 
them into a more complete theory of drag reduction by 
microbubbles. It is the proposition of this Letter how- 
ever that the limit C(z) <§; 1 and d — ► is a relevant 
limit where the theory simplifies considerably and where 
experiments, and especially numerical simulations, can 
give valuable support for the present theory. It is im- 
portant to exhaust the linear effects of drag reduction by 
minute microbubbles before landing on the much more 
involved nonlinear theory. 



tion x is a property of the reference (Newtonian) flow, we 
can take it from Newtonian data. As an example of such 
a calculation we have considered the results of numer- 
ical simulations for a Newtonian channel flow available 
in |2Jj, where the profile S + (z + ) is provided. We have 
used it to compute the response function x and its two 
contributions according to Eqs. I|21 |) -Q22 |) . The results 
are presented in Fig. |2] We see that the qualitative pre- 
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